function nq = cal_nq(const, q_indicator)

    % calculate n(q)
    nq = const.density_omega'*q_indicator;

end